sigma1 = [2,1; 1,2];
sigma2 = [1,0; 0,1];
sigma3 = [1.5,0; 0,1.5];
mu1 = [0;3];
mu2 = [0;0];
mu3 = [5;0];

X1 = gauss_sample(mu1, sigma1, 100)';
X2 = gauss_sample(mu2, sigma2, 100)';
X3 = gauss_sample(mu3, sigma3, 100)';
Xcombined = [X1; X2; X3];
yall = [ones(100,1); 2 * ones(100,1); 3 * ones(100,1)];
perturbation = gauss_sample([0;0], [0.1,0; 0,0.1], 300)';
Xall = cell(2, 1);
Xall{1} = Xcombined;
Xall{2} = Xcombined + perturbation;

num_components = 2;
gmmcommittee = gmmcommittee_train(Xall, yall, num_components);
yall2 = gmmcommittee_classify(Xall, gmmcommittee);
ycorrect = yall == yall2;

figure
hold on
scatter(X1(:,1),X1(:,2),'x','b');
scatter(X2(:,1),X2(:,2),'x','k');
scatter(X3(:,1),X3(:,2),'x','c');
plot_gauss_2d(mu1, sigma1, 1, 'b');
plot_gauss_2d(mu2, sigma2, 1, 'k');
plot_gauss_2d(mu3, sigma3, 1, 'c');

figure
hold on
scatter(Xall{1}(:,1),Xall{1}(:,2), 'x', 'g');
scatter(Xall{2}(:,1),Xall{2}(:,2), 'x', 'r');

figure
hold on
scatter(Xcombined(ycorrect,1),Xcombined(ycorrect,2),'x','g');
scatter(Xcombined(ycorrect == 0,1),Xcombined(ycorrect == 0,2),'x','r');
plot_gauss_2d(mu1, sigma1, 1, 'b');
plot_gauss_2d(mu2, sigma2, 1, 'k');
plot_gauss_2d(mu3, sigma3, 1, 'c');

figure
hold on
scatter(Xcombined(ycorrect,1),Xcombined(ycorrect,2),'x','g');
scatter(Xcombined(ycorrect == 0,1),Xcombined(ycorrect == 0,2),'x','r');
for i = 1 : 2
    for k = 1 : num_components
        plot_gauss_2d(gmmcommittee.classifier{i}.model{1}.mu{k}, gmmcommittee.classifier{i}.model{1}.sigma{k}, 1, 'b');
        plot_gauss_2d(gmmcommittee.classifier{i}.model{2}.mu{k}, gmmcommittee.classifier{i}.model{2}.sigma{k}, 1, 'k');
        plot_gauss_2d(gmmcommittee.classifier{i}.model{3}.mu{k}, gmmcommittee.classifier{i}.model{3}.sigma{k}, 1, 'c');
    end
end